#delimit;
clear all;
set more off;
set logtype text;
capture log close regressions_collapsed_OLS;
local fileloc = "~/KMS_REPLICATION";

log using `fileloc'/log_files/regressions_collapsed_OLS.txt, replace name(regressions_collapsed_OLS);

** Run regression prep;
do `fileloc'/all_regression_prep_hazard_collapsed.do `fileloc' _all;

*XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX;
*XXXXXXXXXXXXXXXXXXXXXXX REGS XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX;
*XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX;

***** VARIABLE LABELS;
label variable died "Infant Death";
label variable male "Male";
label variable black "Black";
label variable asian "Asian";
label variable hisp "Hispanic";
label variable other_race "Other Race";
label variable HS_grad "Mother is HS Grad";
label variable college_grad  "Mother is College Grad";
label variable twins  "Twins";
label variable trip_or_more "Triplets or More";
label variable age19_25 "Mother Age 19 - 25";
label variable age26_30 "Mother Age 26 - 30";
label variable age31_35 "Mother Age 31 - 35";
label variable age36up "Mother Age > 35";
label variable medicaid "Medicaid";
label variable care_first_tri "Care 1st Trimester";
label variable low_weight "Low Birth Weight";
label variable premature "Premature";
label variable weekly_pm10 "PM10";
label variable weekly_co "CO";
label variable weekly_oz "O3";

*XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX;
*XXXXXXXXXXXXXXXXXXXXXXX REGS XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX;
*XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX;

cap program drop effects;
program define effects;
	estadd scalar deaths_compare_`2' = (_b[`1']/1000) * 52 * 100000;
	estadd scalar deaths_compare_`2'_w = ((_b[`1']/1000)*$`1'_within_sd) * 52 * 100000;
	estadd scalar deaths_compare_`2'_b = ((_b[`1']/1000)*$`1'_between_sd) * 52 * 100000;
end;

** Get standard deviations;
preserve;
bysort mother_zip merging_week: keep if _n == 1;
xtset mother_zip;
foreach var in weekly_co weekly_pm10 {;
	xtsum `var';
	global `var'_within_sd = r(sd_w);
	di $`var'_within_sd;
	global `var'_between_sd = r(sd_b);
	di $`var'_between_sd;
};
restore;

** Make coefficients easier to read;
replace died = died * 1000;

*** KEEP ONLY THOSE WITH AT LEAST ONE TRAFFIC SENSOR IN 15 MILES;
egen maxtotflowbase = max(tot_flow_base), by(mother_zip);
drop if maxtotflowbase == 0;

** Use zip by relevant month fixed effects;
eststo clear;

egen zip_rmonth = group(mother_zip event_month);

egen month_by_year = group(event_month event_year);

xi i.birth_quarter i.month_by_year;

areg died weekly_pm10 weekly_co 

	_I*
	
	c.max_temp
	c.humidS
	c.windspeed
	c.rain
	c.days_with_rain
	c.days_with_fog 

	c.max_temp#c.max_temp
	c.humidS#c.humidS
	c.windspeed#c.windspeed
	c.rain#c.rain

	c.max_temp#c.max_temp#c.max_temp
	c.humidS#c.humidS#c.humidS
	c.windspeed#c.windspeed#c.windspeed
	c.rain#c.rain#c.rain
	
	$triweekly_co $triweekly_pm10

	$controls 
	$agespline [fw = weight], absorb(zip_rmonth) cluster(mother_zip);


effects "weekly_co" "CO";
effects "weekly_pm10" "PM10";

estout using `fileloc'/regs/OLS.xls, 
	replace
	lz
	title(KMS OLS)		
	mlabels(, title)
	label
	cells(b(fmt(%9.4f) star) se(fmt(%9.4f))) 
	starlevels(\text{$^{\ast}$} 0.10 \text{$^{\ast\ast}$} 0.05 \text{$^{\ast\ast\ast}$} 0.01)
	stats(N deaths_compare_CO  deaths_compare_PM10 deaths_compare_CO_w deaths_compare_PM10_w deaths_compare_CO_b deaths_compare_PM10_b, label(Obs. "CO" "PM10"  "CO" "PM10" "CO" "PM10") fmt(%9.0fc %9.2f %9.2f %9.2f %9.2f %9.2f %9.2f));
	

log close regressions_collapsed_OLS;
